#include "bc_weir_bottomElevation.hpp"
#include "bc_weir_funcs.hpp"
#include "field.hpp"
#include "boundaryPatch.hpp"
#include "exceptions.hpp"
#include "formatHelpers.hpp"
#include "unifiedMathFunctions.hpp"
#include "geometry.hpp"

namespace hms
{

#define CHEAD(RETURN_TYPE) \
template<typename MeshType> \
RETURN_TYPE Weir_bottomElevation<MeshType>::

CHEAD() Weir_bottomElevation(
	Field<MeshType>& z,
	const BoundaryPatch& p,
	const Field<MeshType>& d,
	const Field<MeshType>& v,
	scalar alpha
) :
	BaseType( Type::weir_bottomElevation, z, p, d, v, alpha )
{}

CHEAD(void) update() const {
	Field<MeshType>& z ( field() );
	const Field<MeshType>& d ( dField() );
	const BoundaryPatch& p { patch() };

	/* this performs coordinate transformation and the computation of dEdge */
	this->updateBuffers(false);

	const auto& dDom { d.domainValues(p) };
	const auto& zDom { z.domainValues(p) };

	z.ghostValues(p) = (dEdge() >= 0.75*dDom).select(
		zDom - 2 * dEdge().square() / (dEdge() + dDom) * ( 1 - cp() ),
		zDom - 8./7 * (
			2/phys::g * ( 4./3 - dDom/dEdge() )
			+ 9./16 * dDom
			- cp() * dEdge().square() / dDom
		)
	);
}

#undef CHEAD

#define INSTANTIATE_FOR_ALL_MESHTYPES(MeshType) \
template class Weir_bottomElevation<MeshType>;
#include "meshInstantiator.hpp"

} // namespace hms
